Phase transition of light on complex quantum networks 
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Recent advances in quantum optics and atomic physics allow for an unprecedented level of control 
over light-matter interactions, which can be exploited to investigate new physical phenomena. In this 
work we are interested in the role played by the topology of quantum networks describing coupled 
optical cavities and local atomic degrees of freedom. In particular, using a mean-field approxima- 
tion, we study the phase diagram of the Jaynes-Cummings-Hubbard model on complex networks 
topologies, and we characterize the transition between a Mott-like phase of localized polaritons and 
a superfiuid phase. We found that, for complex topologies, the phase diagram is non-trivial and 
well defined in the thermodynamic limit only if the hopping coefficient scales like the inverse of 
the maximal eigenvalue of the adjacency matrix of the network. Furthermore we provide numerical 
evidences that, for some complex network topologies, this scaling implies an asymptotically vanish- 
ing hopping coefficient in the limit of large network sizes. The latter result suggests the interesting 
possibility of observing quantum phase transitions of light on complex quantum networks even with 
very small couplings between the optical cavities. 
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I. INTRODUCTION 

Quantum optics and atomic physics have reached a 
level of control over light-matter interactions which not 
only makes it feasible the emulation of condensed matter 
models [l|, Q , but also inspires new architectures envis- 
aging a future quantum Internet with a desired topol- 
ogy ]!, H| . Hence the potential advantage coming from 
a combined optical and atomic approach is twofold: on 
one hand, being able to control a quantum system that 
simulates another one is a way to realize a special pur- 
pose quantum computer; on the other hand, the possibil- 
ity of manipulating new degrees of freedom (not acces- 
sible in condensed matter systems) motivates the exper- 
imental and theoretical study of new quantum systems, 
with the possibility of discovering new physical phenom- 
ena. In this respect an important outcome, coming from 
the combined experimental investigations of atomic and 
optical systems, is the realization of coupled cavity ar- 
rays interacting with local atomic degrees of freedom 
p|. These systems allow for the controlled interaction 
between trapped atoms and local cavity photons, and 
moreover photons are free to hop between coupled cav- 
ities. Changing the details of the physical setup differ- 
ent many-body models can be realized 0, Q , and one in 
particular is of interest in the present work: the Jaynes- 
Cummings-Hubbard (JCH) model (6I-H2I]. Part of the 
interest in these systems is motivated by the possibility 
to investigate new quantum critical phenomena (l3| , like 
the quantum phase transition of light between a Mott- 
like regime and a superfiuid phase [6H12|. Also of inter- 
est is the possibility to generate quantum simulators that 
naturally access non-equilibrium physics [3] . 

In this work we are interested in an additional degree of 



freedom which optical arrays can provide, i.e. the topol- 
ogy of the network underlying the quantum dynamics. 
While regular lattice structures with short-range interac- 
tions are the typical framework in standard quantum em- 
ulator architectures, fiber-coupled cavities may allow for 
the realization of quantum networks with distant effective 
interactions between local degrees of freedom [H, UH, . 
These effective long-range interactions are an important 
ingredient for the construction of quantum networks with 
complex topology, which is a recent topic of interest in 
the quan tum information and complex network commu- 
nity |r7H3l| . Along this line of research, in this work we 
study the effect of the array topology on the phase dia- 
gram of the JCH model. Motivations come not only from 
the possible experimental realization of such systems, but 
also from a number of results, especially in the classical 
context, underlying the importance of networks' topolo- 
gies. Indeed for classical systems it is well known that 
the topology of the network can significantly change the 
phase diag ram of some models, and their critical behav- 
ior [32i . 33] . On the quantum side previous results on the 
Bose-Hubbard model (3~H - [38| are particularly inspiring 
for the present work. In fact, for this model, it has been 
shown by mean-field arguments that the phase diagram 
depends of the maximal eigenvalue A of the hopping ma- 
trix describing the topology of the network |29|j38j . This 
result is valid both in presence of disorder [38[ and in 
presence of a complex topology [29]. Interestingly in a 
complex random network topology, or in an Apollonian 
network [39l.l40|. the maximal eigenvalue of the adjacency 
matrix diverges with the network size. In [29| it has been 
shown, for the Bose-Hubbard model, that this divergence 
implies a non-trivial scaling of the hopping coefficient, in 
order to not suppress the Mott-insulator phase in the 



thermodynamic limit. 

Here we consider instead the properties of the JCH 
model on complex quantum network topologies. Using 
mean-field theory we characterize the phase diagram of 
the model at T = 0, which presents a phase transition 
between a Mott-like regime and a superfluid phase. We 
demonstrate analytically, and confirm numerically, that 
the phase diagram is non-trivial and well defined only 
if the hopping coefficient k scales as the inverse of the 
maximal eigenvalue A of the hopping matrix, i.e. k cx 4 . 
Furthermore we characterize the scaling of the maxi- 
mal eigenvalue for a number of well known complex net- 
work topologies, showing that in many cases the maximal 
eigenvalue A diverges with the network size N. Therefore 
our results are of general interest for a number of different 
complex topologies, and they imply that for complex net- 
work arrays interesting quantum critical behaviors can be 
observed even with very small couplings between differ- 
ent cavities. 

The rest of the paper is organized as follows: in 
Sec. [H] we introduce the Jaynes-Cummings Hamiltonian 
and some of its properties; in Sec. lIIII we characterize the 
Jaynes-Cummings-Hubbard model, its mean-field solu- 
tion, and we consider the scaling behavior of the hopping 
coefficient for different network topologies; Sec. [IV] is de- 
voted to discussions and conclusions. In the Appendix a 
detailed derivation of the mean field solution is provided. 

II. ATOM-PHOTON INTERACTION IN A 
SINGLE CAVITY 

The standard model describing the interaction between 
a two-level atom and quantized electromagnetic modes 
is provided by the Jaynes-Cummings Hamiltonian. In 
the rotating wave approximation, and assuming a single 
cavity mode, the Hamiltonian is given by 

H JC = eo+a- + ucJ a + P(a + a + <j- a*) (1) 

where e is the atomic transition frequency, tu is the field 
frequency, and j3 is the atom-cavity coupling constant; 
a and are the bosonic lowering and raising opera- 
tors, while a are the atomic lowering and raising op- 
erators of the two level system. The eigenstate of this 
Hamiltonian are polaritons, or dressed states, given by 
a combination of atom and field states. In the base 
{|0,4) , |0,f) , 1 , |l,f) • • •}, the atom states are rep- 
resented in the basis of the eigenstates |4) > It) °f the 
Pauli <J Z operator, while the field states are denoted 
with the number operator's eigenstates \n). The Jaynes- 
Cummings Hamiltonian eigenstates arc given by 

| n, — ) = cos#„ \n,l) — sm8 n \n — 

\n, +) = sm6 n | n,4.) + cos6>„ \n - l,t) (2) 

for every n > 1 , where the angle n is expressed in terms 
of the detuning parameter A = e — u> and is given by 

On = ^ arctan ^ 2 ^" ^ . (3) 
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FIG. 1: The energy non- linear dependence of the spectrum 
E(n, ±) — nu) as a function of the detuning A for /3 = 1. The 
non-linearity effects are stronger for low values of n. 

The eigenvalues associated to these eigenstates are given 
by 



E(n,±) =um+ — ± ] J n f3 2 + — (4) 

In addition to the above dressed state, another eigenstate 
of the system is |0) = |0, 4), when n — 0, with the as- 
sociated eigenvalue Eq = 0. The fundamental state of 
the system should be determined for every fixed value of 
the parameters. If we proceed in this calculation we can 
observe first of all that the ground state will be either 
the state with zero polations \n = 0) and energy Eo = 0, 
or one of the states \n—) associated to the eigenvalues 
E n ,-- Indeed, for every fixed number of polaritons n > 1 
we have E n .+ > £"„.-. If we consider the spectrum in the 
limit uj 3> |A|,/3 the state with zero polaritons |n = 0) 
will be the ground state. As we decrease ui we will find a 
point in the parameter space where Eq — , precisely 
at w = A/2 - ^(A/2) 2 + /3 2 . Lowering the value of w 
further we will find a full set of degeneracy points given 
by (for n > 1) 

The energy spectrum of the atom-cavity system, given 
by Eq. |4]), has a nonlinear dependence on n (see the 
energy spectrum in Figure [1]). This anharmonicity in the 
splitting of the energy eigenstates gives rise to nonlinear 
phenomena at the single-photon level. One of the most 
relevant of these is photon blockade, where the presence 
of one photon stops further absorption of photons from 
a coherent light source (4lM43| . 

III. JAYNES-CUMMINGS-HUBBARD MODEL 

Optical cavities, with trapped atoms, can be arranged 
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in arrays where the overlap between different cavities 
wave-functions allow photons to hop from one site to 
another. The Hamiltonian describing this new physi- 
cal scenario is now known as Jaynes-Cummings-Hubbard 
model. The inclusion of optical fibres, or other optical 
devices, can be used to realize more complex geometries, 
where the hopping is not restricted to nearest neighbour 
cavities on a regular lattice d, UH, [l6| . This is precisely 
the kind of situation that we want to investigate in this 
paper. 

To tune the number of polaritons in each cavity a 
chemical potential \i might be used, hence the full JCH 
model will be described by the following Hamiltonian, 

h jch = J2 [Hf c - fiNi] + H h °P (6) 

i 

where 

H- C = ea+ar + cojoi + + <J~a\) (7) 

is the Jaynes-Cummings Hamiltonian for a single cav- 
ity, and Ni indicates the number of polaritons in each 
cavity (Ni = cr^ ' + aja,); while H hop is the hopping 
term. Note that the chemical potential \i is not an ex- 
perimentally tunable parameter for this system. In real 
experiments appropriate preparation schemes have to be 
devised in order to obtain states with different polariton 
number. Recently, using the Rabi model, it has been 
shown that the inclusion of counter rotating terms can 
stabilize finite-density quantum phases of correlated pho- 
tons without the use of a chemical potential [j| . The last 
hopping term in Eq. ((B]) is characterized by a the strength 
K and the adjacency matrix of the underlying quantum 
network t, and is given by 

H hop = -K^T i j(a]a J +a]a i ). (8) 

Let us consider first two extreme cases: one in which 
the hopping strength is very small, and the other where 
the atom-photon interaction is negligible. In the atomic 
limit k//3 <C 1, the Hamiltonian H JCH becomes, to first 
order approximation, the sum of single cavity Hamilto- 
nians H = Hi , with Hi given by 

Hi = H( c - liN. (9) 

The eigenstates of the single cavity Hamiltonian are given 
by Eqs. © for n > 1 , and |0,|) = |0) for n = 0. The 
corresponding eigenvalues are 

££ ± = (u-f i )n+j ±JnP + ^, (10) 

for n > 1, and Eq = for n = 0. The ground state of 
the system can be calculated similarly to the case of a 
single cavity. Indeed, for every cavity we will found that 
the ground state is constituted either by the eigenstatc 
\n, — ) or by the eigenstate |0). 



In the hopping dominated limit k/ [3 3> 1 we can treat 
perturbatively the atom-photon interaction. H JCH re- 
duces, to first-order approximation, to a tight-binding 
hamiltonian H tb given by 

H tb = - fJ,)a\ai - «i T ij( a l a J + a ] a i)- ( n ) 

i " i,j 

The eigenvalues of H tb depends in a simple way from the 
eigenvalues A„ of the adjacency matrix of the quantum 
network: 



E n = N(uj -n- K\ n ). (12) 

The above equation reveals an instability of the system 
for 

kA > uj - [i (13) 

where A is the maximal eigenvalue of the adjacency ma- 
trix t. From this result we can already conclude that the 
maximal eigenvalue of the adjacency matrix set an im- 
portant scale for the strength of the hopping coefficient 

K. 



A. Mean-field theory 

In order to explore the phase diagram of the Jaynes- 
Cummings-Hubbard model we make use of the mean- 
field treatment of the hopping term, which reduces to 
the following approximation 

(Ha) ~ (a,i)a) + a, (a]) - (a,) (a]). (14) 
Therefore the hopping term becomes 

H hop = - K Tl 3 + ai ^3 +K Y1 Tij^J ' ( 15 ) 
id i-.j 

where we have indicated by tpi the local order parameter 
ipi = (a,) (also equal to (aj), due to the gauge symme- 
try of the model). This Hamiltonian displays a phase 
transition between a Mott-Insulator phase, where ipi = 
Vi, and a superfhud phase. In order to study the phase 
diagram of this model, within the mean-field approxima- 
tion we treat H hop as a perturbation and we calculate 
ipi self-consistently, to first order in k, obtaining (see the 
Appendix for more details) 

Ipi = K-^TjjlpjRn, (16) 
3 
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with R n given by 



Rn>\ 



cos 9i 
~E?7 



sin 2 6>i 
^1+ 



(17) 



| \Jn + 1 cos 9 n cos 8 n+ i + y/n sin 9 n sin n+ i | 



n- 

\\/n + lcos# n sin0 n +i — y / risin( 



, COS I 



r n+l| 



^("+1)4- 



+ 



\\/ncos9 n cos0„_i + — 1 sin 9 n sin# n _i| 



E? u 
(n-1)- 



|-\/ncos#« sin(9„_i — i/n — 1 sin cos(9„_i| 

+ _ TP A 1 

^n- £/ (ri-l) + 

where the integer n > depends on the systems pa- 
rameters and it is chosen minimizing the on site energy 
E£_ given by Eq. (fTU|) . If we diagonalize Eq. (fTo]) . along 
with the eigenvalues of the adjacency matrix r, we get 
that the critical line for the transition between the Mott- 
insulator phase and the supcrfluid phase is given by 



nAR n = 1, 



(18) 



where A is the maximal eigenvalue of the adjacency ma- 
trix r. This clearly shows that the phase diagram of the 
model depends on the product kA. On regular graphs 
we have ipi = tp Vi, and A = z, where z is the con- 
nectivity of the lattice. On the other hand for complex 
topology the maximal eigenvalue of the hopping matrix 
A can be significantly different from the average connec- 
tivity of the networks. In particular, for a large variety of 
networks the maximal eigenvalue of the adjacency matrix 
A diverges with the network size. This suggests that, in 
order to have a non-trivial phase diagram for the Jaynes- 
Cummings-Hubbard model, the hopping strength k must 
scale as 
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(19) 



In the following section we will investigate for different 
network topologies the respective scaling, with the net- 
work size, of the maximal eigenvalue of the adjacency 
matrix [4J, |45|. We note here that at the mean- filed 
level the phase boundary is given by Eq. ([15} . With re- 
spect to the phase diagram on regular lattices, the effect 
of the complex topology is the substitution of the aver- 
age degree z of the lattice with the maximal eigenvalue 
A of the network. Therefore the dependence of the phase 
boundary on the detuning parameter A is similar to the 
one observed in regular lattices In fact, as soon as 
the detuning is different from zero, the Mott lobes with 
mean polariton number greater than one are reduced in 
size and shifted to smaller value of the chemical poten- 
tial. This effect is independent on the sign of the detun- 
ing parameter. We remark that on complex network as 
in regular lattices the thermal fluctuations destroy the 



Mott insulator phase. Therefore at finite temperature 
the phase diagram should be composed by a supcrfluid 
regime and a normal fluid. 



B. Regular networks 

For regular networks and regular lattices with connec- 
tivity z the maximal eigenvalue of the adjacency matrix 
A is independent on the network size. In particular, we 
have 



A 



(20) 



In this case the critical line Eq. (|T8)) coincides with the 
one found in the literature using the mean-field approxi- 
mation 0-11]. 



C. Random graphs 

For random Erdos-Rcnyi graphs with finite connectiv- 
ity and Poisson degree distribution, it has been proven 
[43 that 



A OC \f~kn 



(21) 



where k max is the maximal degree of the system. For 
random networks with a finite connectivity we have 



In N/ In In N, therefore 



A(N) cx 



In AT 
In In N ' 



(22) 



Considering the scaling given by Eq. (fT!?]) . we have that 
the hopping strength has to satisfy the following relation 
in order to have a non-trivial phase diagram 



k(N) cx 



In In N 
In N ' 



(23) 



D. Random scale-free networks 

Scale-free networks provide one of the most interest- 
ing and most studied topology for the analysis of phase 
transitions occurring on them. Indeed, classically, on 
scale-free networks with degree distribution P(k) cx A: -7 
and power-law exponent 7 € (2, 3] the phase diagram of 
the Ising model |4q-|49{. and the percolation transition 
[50l . HU change drastically due to the diverging second 
moment of the average degree (k 2 ) . A similar observation 
can be also made for the epidemic spreading model on 
annealed complex networks (52[, i.e. complex networks 
in which the links are dynamically rewired. Moreover, 
the spectral properties of the complex networks deter- 
mine the critical behavior of the epidemic spreading on 
complex quenched networks [H, [54 1 for the 0(N) model 
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FIG. 2: (Color online) Mean-field phase diagram of the 
Jaynes-Cummings-Hubbard model with A = on a random 
scale-free network with power-law exponent 7 = 2.2 for dif- 
ferent value of N, i.e. N = 100 (top panel) N = 1000 (middle 
panel) TV = 10000 (bottom panel). The phase diagram scales 
with the maximal eigenvalue that is given by A = 5.98 (top 
panel) A = 11.07 (middle panel) and by A = 23.42 (bottom 
panel). The solid lines denote the analytic perturbative solu- 
tion in mean field of the model. 





w 



FIG. 3: (Color online) Left panel: the first 3 generations of 
the Apollonian graph. The nodes and links added to construct 
the 1st, 2nd and 3rd generations are shown. Right panel: the 
5th generation Apollonian network. 
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and the stability of the synchronization dynamics 



In random scale-free networks, with power-law degree 
distribution P(k) oc k 1 it has been proven [44[ that the 
maximal eigenvalue scales like 



A 



oc 



y/kmax for 7 > 2.5 
W 7 < 2.5. 



Moreover, the maximal degree of the network satisfy 
kmax = min [N 1 ^ 2 , TV 1 /*- 7-1 )] , where we have considered 
the structural cutoff of the degrees of the network for 
A < 3. Therefore, the maximal eigenvalue of the net- 
work A follow a different scaling with the network size, 
depending of the power-law exponent 7, 

f #l/[2(7-i)] fo r 7 > 3 

A(N) oc I N 1 / 4 for 2.5 < 7 < 3 
[ Af(3-7)/2 for 7 < 2.5. 

Finally, the hopping coefficient k that ensure a non-trivial 
phase diagram [see Eq. (fl9|) ] scales with the network size 
N according to the following rules 



k(N) 



N -l/[2(y-l)] 
jy-(3- 7 )/2 



for 7 > 3 
for 2.5 < 7 < 3 
for 7 < 2.5 



Figure[2] shows both the analytic perturbative solution in 
mean-field of the JCH model, and the numerical non- 
perturbative mean-field evaluation of the phase diagram. 
As can be seen there is substantial agreement between the 
two, and furthermore one can observe the dependence on 
the size of the network of the detailed location of the 
critical lines. 



E. Apollonian networks 

We consider here Apollonian networks [39| which are 
constructed through a 2D Apollonian packing model in 
which the space between three tangent circles placed on 
the vertices of an equilateral triangle is filled by a max- 
imal circle. The space-filling procedure is repeated for 
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every space bounded by three of the previously drawn 
tangent circles. The corresponding Apollonian network is 
constructed by connecting the centers of all the touching 
circles [Figure [3] (left)] . The resulting network is scale- 
free with power-law degree distribution p(k) = Nk~ x 
and A = l + ln(3)/ln(2) ~ 2.585. Also these networks are 
known to have di verg ing maximal eigenvalue A of their 
adjacency matrix [40(. In Figure 0] we plot the maximal 
eigenvalue A of the apollonian network as a function of 
the network size N . We can fit the numerical results with 
the function 



eigenvalue of this network, for p = 1 will increases with 



A(N) oc N' 



0.23 



(24) 



Therefore the hopping coefficient k [that needs to scale 
according Eq. (fT9|) ] scales for large N as 



k(N) oc N 



-0.23 



(25) 




100 r 



A=3.7N 



< 10 r 



FIG. 5: The maximal eigenvalue of the Small- World network 
as a function of p for different network sizes N. In the limit 
p — > the small world model is a regular chain with average 
connectivity z = 4, therefore in this limit A = 4. The data 
are averaged over 100 network realizations. 



N, as in the random graph case [Ea.(|22|l]. while for the 
case p = it will be independent on N, as in regular net- 
works. In Figure [5] we see how the maximal eigenvalue 
of the network changes with N for intermediate values of 
the probability p when the network is small-world. 



\o 100 1000 10000 

N 

FIG. 4: Scaling of the maximal eigenvalue A of the Apollonian 
network as a function of the size N of the network. 



F. Small- world networks 

Small- world networks structures characterize a system 
that has at the same time a small diameter, like ran- 
dom networks, but also has a high clustering coefficient, 
similarly to regular networks [591 ]. In particular, we can 
follow the construction proposed in [59( : we start from a 
regular chain in which each node is linked to the nearest 
neighbours and to the next-nearest neighbours; then ev- 
ery link is rewired with probability p to another random 
node of the network. For p = the small-world net- 
work is a regular lattice in one dimension; for p = 1 the 
small-world network becomes one instance of a random 
graph network; finally, for every intermediate value of p 
we observe the small-world network with small average 
diameter and a high clustering coefficient. The maximal 



IV. CONCLUSIONS 

In conclusion we have studied the Jaynes-Cummings- 
Hubbard (JCH) model on complex quantum networks. 
We have shown that the phase diagram derived in the 
mean-field approximation, depends crucially on the max- 
imal eigenvalue A of the hopping matrix. In particular 
the phase diagram depends on the product kA. This 
implies that in order to have a well defined phase dia- 
gram in the large network limit, the hopping coefficient 
k should scales proportionally to 1/A. The eigenvalue A 
is equal to the connectivity of the network, for regular 
networks and lattices, but for complex random networks 
it generally increases with the network size. In this paper 
we have listed for a large class of networks the scaling of 
the maximal eigenvalue A with the network size N. For 
complex networks that have a diverging A, the hopping 
coefficient should be a decreasing function of N in order 
to observe the phase transition from the Mott-like regime 
to the supcrfluid phase. This result implies the possibil- 
ity of observing quantum critical behaviours in arrays 
whose cavities are weakly coupled, assuming the effec- 
tive realization of the proper complex quantum network 
topology. 
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Appendix A: Mean field solution of the JCH model 

The Jaynes-Cummings Hamiltonian (we set h = 1) 

H JC = 6<t + <j- +uja^ a + I3(a + a + 0--^) (Al) 

is obtained in the rotating wave approximation and in 
the limit /3 -C e, oj. The total number of excitations is a 
conserved quantity, and it is given by the sum of electro- 
magnetic and atomic excitations N = a^a + (j + a~ . The 
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interacting part of the Hamiltonian connects n-sectors 
which differs only by one photon excitation 



1) It) 



\n) ||) 



(A2) 



The JC Hamiltonian can then be block-diagonalized in 
different sectors, each labelled by n, and each sector 
spanned by {\n — 1) |t) , |") Choosing this set as 

the basis for the n— th sector, Eq.([2|) and Eq.(j4|) in the 
main text provide the expressions for the eigenvectors 
and eigenvalues of the JC Hamiltonian. 

Considering the situation of a network of cavities, 
whose coupling is effectively described by an hopping 
term, we have the Jaynes-Cummings-Hubbard Hamilto- 
nian described in Eq.([6]). As explained in the main text, 
the mean-field treatment of the hopping term allows us 
to approximate the JCH Hamiltonian as follows 

h mf = J- ea+ar + ujaja, + P{o+ ai + a~aj) 

i 



i -J 



•■J 



Considering the hopping term as a perturbation to the 
atomic limit, the order parameter of the model is pro- 
vided by ipi = (ai) gs , where the expectation value is cal- 
culated with respect to the ground-state of the Jayncs- 
Cummings- Hubbard model to first order in perturbation 
theory. Note that the order parameter can be assumed 
real due to the gauge symmetries of the Hamiltonian Q . 
The self-consistent equation for the order parameter can 
then be written as 



ipi = (ni\ai\ni), 



(A3) 



where |?ii) = |n, — )° + In) 1 is the approximation of 
the ground-state to first-order in the perturbation, while 
\n, — )° is the ground-state of the unperturbed Hamilto- 
nian (see Eq[2]in the main text), and 



> 1= E 



\k,a\H 



MF\ 
hop I 



Ff - F/ 



\k,a}°. 



From Eq. (|A3|) we have 



(nx\ai\m) = °(n, -|a,|n, -)° + 1 (n\a l \n) 1 

+ °(n,-\a i \n} 1 + 1 (n^n, -)°. 

Keeping only non-zero terms to first order in k we are 
left only with the last two terms in the above equation. 
First we explicitly calculate 



{n,-\ai\n) = ^ — ^ 



MFl _\0 
hop K°) / o 



k,a=± 



E. 



n,—\a,i\k, a) c 



It is easy to check that the only non-zero terms in the 
sum are given by k = (n + 1), and a = +, — . It follows 
that the non-zero contribution to the expectation value 
of Hffjp is provided only by — T ij^[ V'j- From the 
explicit form for the ground-states of the unperturbed 
Hamiltonian (see Eqj2]in the main text) we obtain 



°(n, -\ai\(n + 1), -)° = VnTTcos6>„ cos6»„ + i 

+y/n sin n sin 6 n +i 



(n, — \ai\(n + 1), +) = y/rT+1 cos n sin n +i 
— yfa sin 6 n cos Q n + 1 . 

We can proceed similarly for the calculation of 
1 (n|cii|n, — )°, obtaining 

((n — 1), — \ai\n, — ) = \fncosO n cosO n -i 

sin n sin 9 n —i 



°((n - 1), +|o,|n, — )° = v^cos^n sin6»„_i 

— \fn — 1 sin n cos n -\. 



Now let us define 

Rn = - E 



\ (n-\ ai \(n+l),an 



a=± 

|°((n-l),a|a 4 |n,-)°| 2 



Putting everything together we have the following self- 
consistent equation to first order in perturbation theory 
(for n > 0) 



Ipi = nR n 2J Tijl/>j 



The case n = can be calculated in the same way and 
one has 



where 



i? = -£ 



l<0,| |o<|l,a)| 5 



a=± ^0,4. a l,a 

with (0,| |a,|l, — ) = cos Oi, and (0,4, |fflt|l, +) = sin#i- 



